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1. INTRODUCTION 

Hysteresis is a phenomenon where the present condition is strongly depended on the past 
circumstance [1]. This trend can be found in functional material such as for example: magnetic material [2], 
[3], ferroelectrics [4] and semiconductors [5]. The hysteresis curves show the response of materials when an 
external stimulation is applied. In magnetics, a hysteresis curve illustrates the behavior of magnetization due 
to the applied external magnetic field. The hysteresis of ferroelectrics expresses the effect of an electric field 
to the polarization. Hence, the study of hysteresis of the functional materials is very important to understand 
the characteristic of the related materials. It is the basic step before functional materials can be utilized in 
many applications such as for example: measurement devices [6], field driven devices [7], LED devices [8], 
etc. 

Magnetoelectric multiferroics are well known as material that possess both ferroelectricity and 
ferromagnetism at the same phase [9]. It means that the materials retain electric polarization and the 
magnetization which determine the properties of material. The ferroelectricity is related to the 
ferromagnetism through magnetoelectric coupling which responsible for the cross behaviour: the application 
of an electric field not only change the electric polarization in ferroelectrics phase, but also disturb the 
magnetization in magnetic phase [10]. In other way, the magnetic field can be used to drive the electric 
polarization [11]. As a consequence, there are four types of hysteresis for magnetoelectric multiferroics 
system. First, it is magnetic hysteresis relating magnetization (M) and magnetic field H. Second type is the 
electric hysteresis connecting the electric polarisation P with the electric field E. The third and fourth are 
magnetoelectric (cross-type) hysteresis linking M with E and P with H. Considering the potential application 
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of multiferroic materials (for example: electric field-controlled magnetic data storage [12]), the study of the 
hysteresis of this type of materials is unavoidable. 

A numerical hysteresis loop for ferroelectrics is obtained by solving Landau-Khalatnikov (LKh) 
equation of motion which describes relaxation processes near phase transition [13], while numerical 
hysteresis for magnetic (ferromagnet) is acquired using Landau-Lifshitz-Gilbert (LLG) equation of motion 
which illustrates precession of the magnetization with the damping is considered [14], [15]. Hence, the 
numerical hysteresis for magnetoelectric multiferroics can be calculated using the electric LKh equation of 
motion and the magnetic LLG equation of motion simultaneously. However, the calculation by combining 
LKh and LLG equations of motion such in [16] is tough since LLG model requires more than one equation. 
Since obtaining hysteresis curve is basic step to understand the characteristic of the materials, it is better if 
the calculation to get hysteresis is less complicated. 

Motivated by the facts that calculation using LKh theory is simpler than calculation employing LLG 
equation of motion and also there are some studies in magnetic system using approximation based on LKh 
equation of motion [17], [18], we propose in this numerical study to approach the problem of obtaining the 
hysteresis of the magnetoelectric multiferroics by employing Landau-Khalatnikov (LKh) theory for both 
systems (the ferroelectric and magnetic). In this paper, we want to show that this treatment is appropriate 
enough to give the initial feature, especially for obtaining the cross hysteresis in magnetoelectric system. As 
long as we know, there is no report regarding this treatment in multiferroics. 


2. RESEARCH METHOD 
The study begins by setting the energy density of the materials. Since magnetoelectric multiferroics 
comprise of ferroelectric and ferromagnet with magnetoelectric coupling, the energy density can be written 


where Fr, Fm and Fme represent the energy density for ferroelectric, ferromagnet phase and magnetoelectric 
interaction. Based on Landau phenomenological theory, the energy density of a system can be written as a 
polynomial of the order parameter of material. Here, the energy density for ferroelectric is represented by 
Ginzburg-Landau equation with electric polarization (P) as order parameter expanded until quartic form as: 


1 1 
Fg = 5aP* rs ie EP (2) 
with a is temperature dependent constant defined as a = ag(T — Tç). Here, ao and b in equation (2) represent 
dielectric stiffness constants while E is an external electric field. Using similar procedure, we determine the 
energy density for ferromagnet with magnetization (M) as the order parameter, in the form 


Fy = aM? + ~BM4 — HM (3) 


where œ = &o(T — T?") is time dependent with œo and £ are magnetic stiffness constants. Here T¥ is Curie 
temperature for magnet and H represents an external magnetic field. According to symmetry argument, the 
energy density represents bi-quadratic coupling such as [19], [20] 


Fug = —5KP?M? (4) 


where x represents magnetoelectric coupling constant, should always exist in all magnetoelectric multi- 
ferroics. The next step is performed by substituting the form of density energy equation (1) into the Landau- 
Khalatnikov equation of motion 


Ot = —y (5) 


where y is phenomenological damping constant and O; represents the order parameter of the system. Since 
magnetoelectric multiferroics have electric polarization (P) and magnetization (M) as the order parameters, 
this step yields two coupled equations such as 


its aP + bP? — kPM? —E (6a) 
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for ferroelectric phase, and 
——— = gM + fM? — kP?M — H (6b) 


for ferromagnet phase. 

Using rescalling processes, we can write LKh equation of motions in equation (6) with the 
dimensionless polarization p and dimensionless magnetization m. Then, the rescalled equations can be 
written in the form which is appropriate for numerical calculations as 


(1 — yaAt)Pi+ı + yaAtpi,ı + yaAtKm}, pis. — Pi — YAACE {41 = 0 (Ta) 
and 
(1 — yma At)Mi41 + Ymadtm3,, + yma AtK p? 1 Mi+ı — M; — yadthj,, = 0 (7b) 


where P=Pp, M=M,m with P, and M, represent the spontaneous electric polarization and lattice 
magnetization. The external electric field, the magnetic field and ME coupling are modified as, 


1 1 
ESE (=) ,m=H (4) and K = (5) k. Here, the subscript of order parameter is related to the time as 


E; = E(t) and E;,, = E(t + At). The last step is solving the coupled LKh equation in equation (7a) and (7b) 
numerically. In this process, the electric (p — E) hysteresis and cross (m — E) hysteresis can be obtained by 
setting the external magnetic field to zero while variating the electric field. On the other hand, by setting the 
external electric field to zero while variating the magnetic field, we can get magnetic (m — h) hysteresis and 
cross (p — h) hysteresis. 

In this research, we also study the role of magnetoelectric coupling in affecting the shape of cross 
hysteresis curve. It can be performed by solving the hysteresis formulas in several values of the 
magnetoelectric coupling. In practice, changing the values of the coupling is very difficult. However, the 
result can be used to predict the strength of the magnetoelectric coupling when the measurement of the cross 
hysteresis is performed. 


3. RESULTS AND ANALYSIS 
In this numerical calculation, we set > = 101° statC?/cm*; - = 10° Oe’, typical for ferroelectric and 


ferromagnetic phase. We also keep the value At = 0.05; ymaæAt = 0.005 and scalled magnetoelectric 
constant K = 5 x 103. Solutions of the LKh equations of motion in equation (7) are illustrated in Figure 1. 
We obtained four types of hysteresis: (p-e), (m-h), (p-h) and (m-e) hysteresis. The hysteresis for (p-e) in 
Figure 1(a) and (m-h) in Figure 1(c) have well known symmetrical shape where the electric polarisation p 
and the magnetic magnetization m have range from negative to positive values. Hence, both types of 
hysteresis posses coercive fields which are required to bring the electric polarization or magnetization to zero 
value. It means that the electric polarization p is able to be flipped into the opposite direction using an 
external electric field e, while the magnetization m can be flipped by applying an external magnetic field h. 
These (p-e) and (m-h) hysteresis curves have similar shape to the hysteresis for ferroelectric materials and 
ferromagnet materials from the earlier researches. 

The calculated cross (m — £) and (p — h) hysteresis which is the most important finding are shown 
in Figure 1(b) and (d). We found that the shape of hysteresis curves is very diferent compare to the shape of 
(p-e) and (m-h) hysteresis. These cross (m — £) and (p — h) hysteresis posses butterfly shape with the range 
of the electric polarization or magnetization is only involving positive values. These numerical results 
qualitatively agree to the results from the previous experimental reports [21], [22]. The butterfly hysteresis 
curves in Figure 1(b) and (d) do not have coercive field since the electric polarization or magnetization could 
not be brought to zero value. It means that the electric polarization and magnetic magnetization in these 
butterfly shape curves in Figure 1(b) and (d) could not be flipped into the opposite direction. We argue that it 
occurs since bi-quadratic type of magnetoelectric coupling does not sensitive to the direction of the applied 
field. Both of the magnetization and the electric polarization are in quadratic term, means that the strength of 
that parameters which is important and not the direction of the electric polarization or magnetization. 
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Figure 1. The hysteresis of magnetoelectric multiferroics with bi-quadratic couplings. Figure (a) and (c) are 
illustrating the electric (p-e) hysteresis and the magnetic (m-h) hysteresis. The other curves showing “cross” 
hysteresis where (b) displays the response of magnetization to the applied electric field while (d) indicates the 
response of polarization to the applied magnetic field 


We also study the role of the strength of magnetoelectric coupling in determining the cross- 
hysteresis in multiferroics materials. The results are shown in Figure 2 with the strength of the coupling in 
Figure 2(a) and Figure 2(c) is weaker than the coupling in Figure 2(b) and Figure 2(d) by factor 107. It can 
be seen that the response of the magnetization or the electric polarization to the applied electric or magnetic 
fields are weak when the strength of the magnetoelectric is small. It is clear when the electric/magnetic field 
is applied to the magnetoelectric (ME) multiferroics, it will change the state of the electric/magnetic system 
significantly. If the ME coupling is small, then these changes of the electric/magnetic system will only give 
the small disturbance to the magnetic/electric phase which in turn resulting weak responses of the 
magnetization/electric polarization. This agrees with the general feature of magnetoelectric materials: the 
cross effect is hardly to sense if the magnetoelectric coupling is weak. 

Overall, the results stated that the approach using Landau-Khalatnikov equation is able to produce a 
qualitatively good prediction. The shape of cross hysteresis evidently similar to the previous experiment 
study. The approximation also yields reasonably feature regarding to the strength of the coupling. 
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Figure 2. The cross hysteresis of magnetoelectric multiferroics with bi-quadratic couplings. Figure (a) and (b) 
are showing the cross (m-e) hysteresis. The bi-quadratic coupling in (a) is smaller by factor 10 than the 
coupling in (b). The similar patern is also illustrated in (c) and (d) for cross (p-h) hysteresis. Here, the 
strength of magnetoelectric coupling in (c) is 10? times weaker than that in (d) 


4. CONCLUSION 

In this research, we prove that Landau-Khalatnikov equation of motion may be used entirely for 
both ferroelectric and ferromagnetic phase to approximate the hysteresis of the magnetoelectric multiferroics 
with bi-quadratic phase. The approximation is able to produce the ‘butterfly-shape’ for cross (p-h) and (m-e) 
hysteresis which agree to the previous experiment result. Hence, this approximation is adequate to initial 
study to get the outline feature of the system before we perform the detail treatment using Landau-Lifshitz- 
Gilbert equation for spin dynamics in magnetic system. 
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